S-wave anisotropy estimate by automated image registration

ABSTRACT

The present disclosure provides a system and method for estimating fracture density within a subsurface formation from S-wave seismic data. In one embodiment, the S-wave seismic data is separated into fast (“S 1 ”) and slow (“S 2 ”) data. A computer is used to compute local similarity of the S 1  and S 2  data and to compute a cumulative time-difference by which the S 2  data lags the S 1  data from the local similarity. Based on the computed cumulative time-difference, the fracture density of a subsurface formation is estimated.

CROSS-REFERENCE TO RELATED APPLICATION

This application is the National Stage of International Application No.PCT/US2012/037583 that published as WO 2013/048508, filed May 11, 2012,which claims the benefit of U.S. Provisional Application No. 61/539,302,filed Sep. 26, 2011, entitled S-WAVE ANISOTROPY ESTIMATE BY AUTOMATEDIMAGE REGISTRATION, each of which is incorporated herein by reference,in its entirety, for all purposes.

FIELD OF INVENTION

This invention generally relates to the field of seismic prospecting andreservoir delineation and, more particularly, to fracturecharacterization in a subsurface through processing of seismic data.More specifically, the invention involves determining a time-continuousS-wave time-difference volume in the subsurface to estimate or predictfracture intensity (or density) within a subsurface formation.

BACKGROUND

Usually, fracture networks, especially in tight-gas sands, are exploitedfor efficient hydrocarbon recovery from the reservoirs [1, 2].Sometimes, hydrocarbon recovery completely relies on the exploitation ofthe natural fracture networks in the subsurface.

Several geophysical techniques are available for characterizing fracturenetworks in the subsurface and each has its own advantages anddisadvantages. All these techniques can be divided into two broadcategories: (1) direct measurements and (2) indirect (or remote)measurements. An example of a direct measurement is a well-bore basedmethod. Usually a geophysical instrument is sent into the well-bore andthe geophysical tool measures the subsurface properties such as seismicvelocities. These subsurface data are used to predict the fractureproperties of the subsurface [3]. Although these types of techniques arevery reliable, they provide fracture properties only at the well-borelocation. Away from the well-bore, these methods cannot be trusted forfracture characterization.

An example of an indirect or remote measurement is surface seismicmethod. Surface seismic methods are one of the most common techniquesfor subsurface imaging. Seismic P- and S-waves are the two types ofseismic waves that are used for this purpose. A P-wave source such asdynamite is used to excite P-wave energy which travels down thesubsurface and reflects back both as P- and S-waves. These reflectedwaves are captured by surface receivers. These reflected energies areused to generate subsurface images and to derive other subsurfaceproperties. P-waves are recorded by vertically oriented receivers andS-wave energies are recoded by horizontally oriented receivers. Thereflected P-wave energies are traditionally called PP modes and thereflected S-wave energies are called PS or converted-wave modes.

In the past, geophysicists have proposed and implemented a number oftechniques to characterize fractures using surface seismic data.Fractured reservoirs are known to behave as an azimuthally anisotropicmedium on the scale of seismic wavelengths [4]. Ruger and Tsvankin [5]showed that PP-reflectivity in fractured reservoirs varies with thefracture azimuth. They also gave analytical expressions forPP-reflectivity which could be used to estimate fracture intensity (ordensity) of the medium. Methods based on this property of the PP-modeare called AVAZ-based methods.

S-waves, also called shear waves, travelling through a fractured mediumsplit into fast (S₁) and slow (S₂) modes. The particle motions of S₁-and S₂-waves are polarized parallel and perpendicular to fracturestrike, respectively. S-waves polarized parallel to fractures (S₁) havea greater velocity than the S-waves polarized perpendicular to fractures(S₂). The difference between the fast and slow S-wave velocities isdirectly proportional to fracture density; i.e. the larger the fracturedensity, the larger the difference between velocities. This phenomenonis called S-wave birefringence [6]. A number of fracturecharacterization methods have been proposed based on this property ofS-wave [20].

Alford [7, 15] proposed a technique for a vertical seismic profile (VSP)geometry that includes rotating, in a synchronic way, source andreceiver geophone by linearly combining the two polarizations. Themethod requires two orthogonal source components and two orthogonalreceiver components. A 2×2 data matrix is formed and the energy in theoff-diagonal terms are minimized by tensor rotation. The angle at whichthe off-diagonal energy is minimized is the azimuth of the fractures inthe subsurface. The main disadvantage of this method is that theestimated fracture properties are only reliable at the VSP location.

Winterstein and Meadows [8] reported that the subsurface rarely has onlyone fractured layer; instead, many fractured layers with varyingfracture orientations are more common They proposed a coarse-layerstripping technique to deal with this problem. The following is the ideabehind their method; first rotate and find the time-difference betweenS₁ and S₂ for the arrivals from the bottom of the first fractured layer,then subtract the one- or two-way time (depending on whether the data isVSP or surface seismic) from the arrivals from the bottom of nextfractured layer and correct for time lag by the first fractured layer.The procedure is repeated for subsequent fractured layers.

Gaiser [9, 14] extended the method of Alford [7, 15] to characterizesubsurface fractures using surface seismic PS data. Unlike the method ofAlford [7, 15], Gaiser's technique uses surface seismic data forfracture characterization. Gaiser's method can also perform coarselayer-stripping in the presence of multiple fractured layers.

Bansal et al. [18] developed a method to perform true-amplitudelayer-stripping using surface seismic PS data. Unlike the method ofGaiser [9, 14], Bansal's method perform a scan of offset and azimuths ofthe surface seismic data to obtain the optimum dataset to perform layerstripping. This method also produces true-amplitude fast and slowS-waves which can be used to perform seismic inversion to in order topredict the lithology of the subsurface.

All the previous fracture characterization methods require the data tobe divided in several time windows before layer-stripping is performedto compute the slow S-wave time-lag and the fracture orientation. Theestimated S-wave time-lag represents the cumulated traveltime differencebetween the fast and slow S-waves in a given time window. In order tolocate the highly fractured zones, length of the analysis time windowsshould be as short as possible such that the estimated S-wave time-lagrepresents the local anisotropy strength and not a cumulative effect.

In the field data, however, the length of time windows depends on thesignal-to-noise ratio. The time-difference between the fast S-wave (S₁)and slow S-wave (S₂) modes are computed by cross-correlating thereflection events in the given window. If the time windows are small andthe signal-to-noise ratio is not sufficient enough, cross-correlationprocess becomes unstable. This limitation on the conventionallayer-stripping usually forces one to choose large time windows foranalysis. This problem becomes more aggravated in land converted-wave(PS) and pure S-wave data which are notoriously noisy. The noiseprecludes identifying highly anisotropic zones from moderatelyanisotropic or isotropic zones.

Thus, there is a need for improvement in this field.

SUMMARY OF THE INVENTION

The present invention provides an improved system and method forfracture characterization in a subsurface formation using geophysicalseismic data.

In one embodiment, the invention is a computer-implemented method forestimating fracture density within a subsurface formation from S-waveseismic data, the method comprising: separating the S-wave seismic datainto fast (“S₁”) and slow (“S₂”) data; using a computer to compute localsimilarity of the S₁ and S₂ data; using a computer to compute acumulative time-difference by which the S₂ data lags the S₁ data fromthe local similarity; and estimating the fracture density from thecomputed cumulative time-difference.

BRIEF DESCRIPTION OF THE DRAWINGS

The patent or application file contains at least one drawing executed incolor. Copies of this patent or patent application publication withcolor drawing(s) will be provided by the U.S. Patent and TrademarkOffice upon request and payment of the necessary fee.

The present invention and its advantages will be better understood byreferring to the following detailed description and the attacheddrawings in which:

FIG. 1 is a flow chart showing basic steps for the fracturecharacterization according to one embodiment of the present disclosure;

FIG. 2A depicts an exemplary similarity scan and warping function for aseismic pair;

FIG. 2B depicts the cumulative S-wave time difference for the seismictrace information depicted in FIG. 2A;

FIG. 3 illustrates S₁ and uncompensated S₂ traces;

FIG. 4 illustrates S₁ and time compensated S₂ traces;

FIG. 5 illustrates a cumulative S-wave time-difference volume based onthe seismic traces shown in FIGS. 3 and 4;

FIG. 6 illustrates an instantaneous S-wave time-difference volume basedon the cumulative S-wave time-difference volume shown in FIG. 5.

The invention will be described in connection with example embodiments.However, to the extent that the following detailed description isspecific to a particular embodiment or a particular use of theinvention, this is intended to be illustrative only, and is not to beconstrued as limiting the scope of the invention. On the contrary, it isintended to cover all alternatives, modifications and equivalents thatmay be included within the scope of the invention, as defined by theappended claims. Persons skilled in the technical field will readilyrecognize that in practical applications of the present inventivemethod, it must be performed on a computer, typically a suitablyprogrammed digital computer.

DESCRIPTION OF THE SELECTED EMBODIMENTS

The disclosed methodology describes a process to generate atime-continuous time-difference volume between the fast S-wave (S₁) andslow S-wave (S₂) modes. The time-difference between S₁ and S₂ modes isan indicator of fracture intensity within a subsurface formation. As aresult, a larger time-difference indicates a larger fracture intensityor density. This new, additional information will help to locate highlyanisotropic or fractured zones: critical information for well planningin unconventional reservoirs. As used throughout this document, atime-continuous volume is understood mean an attribute volume in whichthere is a defined value at each time sample.

The flowchart of FIG. 1 will be referred to in describing one embodimentof the present disclosure. The depicted process (100) first requiresacquisition of multi-component, multi-azimuth data (101). Azimuth isdefined for a particular source-receiver combination. The direction(relative to true North or some other reference direction) of the lineconnecting a source-receiver pair is called the azimuth of thatparticular source-receiver pair and associated seismic data.Traditionally, only the vertical component of the seismic wavefield,which is dominated by the P-wave energy, is acquired. For certainapplications, all three vector components of the wavefield are alsoacquired (using a motion-detector type of seismic receiver).

In this type of acquisition, a P-wave source is used which may be, butis not limited to, dynamite, a vertical vibrator, or air gun. Thevertical component of the data mostly contains P-wave energy and the twohorizontal components carry converted-wave PS energy. The PS energy isdefined as the P-wave energy reflected back from a reflector as S-waveenergy; i.e., P-wave goes down and some of that energy is reflect backup in an S-wave mode.

The process continues as the acquired PS energy is rotated or resolvedinto radial (direction along the line connecting source to receiver) andtransverse (direction perpendicular to the line connecting source toreceiver) components (103). Free-surface related seismic noise, such assurface-waves and free-surface multiples, is then removed from theradial and transverse components (105). After noise correction, normalmoveout (NMO) correction is applied (107) on the data to flatten thereflections. In other embodiments, the reflections may be flattened bypre-stack time migration. Though not depicted in FIG. 1, someembodiments may include the additional step of stacking the flatteneddata. By stacking the data, the signal-to-noise ratio may be improved.As appreciated by those skilled in the relevant art, steps 101-107 arestandard processing steps and are routinely applied in seismic dataprocessing.

The two horizontal components of the recorded wavefields, whichprimarily have S-wave reflected energies, carry a mixture of fast andslow S-waves. The S-wave data recorded over a fractured subsurface isseldom separated in pure S₁ and S₂ modes. In order to do so, thehorizontal components (radial and transverse components) first need tobe rotated to fast and slow S-wave directions before S₁ and S₂ imagescan be time-aligned. The way this is done is dependent upon whetherthere is a single dominant fracture direction in the subsurface. Ifthere is only a single dominate fracture direction in the subsurface(109), the horizontal components can be simply vector-rotated tofracture strike and normal directions to produce S₁ and uncompensated S₂seismic sections (111).

However, if the fracture orientation varies with depth (113), it is notpossible to apply a single vector rotation on the whole section. In thiscase, layer-stripping is performed to produce an S₁ seismic section anda registered (or time-compensated) S₂ section (115). A number of methodshave been published on layer stripping from surface seismic data. Toname a few, Gaiser [9, 14] published a method called “3-D convertedshear wave rotation with layer stripping”. Another method was publishedby Thomsen et al., [11] called “coarse layer stripping of verticallyvariable azimuthal anisotropy from shear-wave data”. Granger et al.,[12] developed a method to find the fast S-wave direction whichcorresponds to the fracture orientation. Haacke et al. [17] proposed amethod of layer-stripping in marine data. Crampin [13] gave a detaileddescription of S-wave propagation in fractured media which led todevelopment of the layer-stripping technique. Bansal et al. [18]developed a method to perform true-amplitude layer-stripping usingsurface seismic PS data. Bansal's method performsf a scan of offset andazimuths of the surface seismic data to obtain the optimum dataset toperform layer stripping.

Conventional layer-stripping techniques typically consist of fourprimary steps: (1) division of data in several layers (or windows), (2)rotation of the radial and transverse components to fracture strike andnormal in first window to generate fast and slow S-waves, (3)cross-correlation of the fast and slow S-waves to estimate thetime-difference between the two in the first window and (4) shifting theslow S-wave by the estimated time-difference to align the fast and slowS-wave events or reflections in first window. This process is performedmultiple times in a top-down fashion revealing fracture orientation andS-wave time-difference of each subsequent fractured layer (or window).This process produces pure S₁- and registered S₂-mode seismic sections,with the registered S₂-mode seismic sections resulting from the timeshifting noted above in step 4. This time shifting is necessary so thatthe fast and slow waves can be accurately determined in the next layer(or window). In one non-limiting embodiment, the method proposed byBansal et al. [18] is used as the layer stripping technique. However, itmust be noted that the embodiments of the present disclosure are notdependent on the type of layer stripping method used.

After the S₁- and S₂-modes have been generated, the process continues bydetermining the time-difference between the two S-waves. Fomel and Jin[19] previously developed a method to register time-lapse seismic imagesusing the local similarity attribute which provides a smooth continuousmeasure of similarity between two images. Their solution produces apoint-by-point aligned time images, time-continuous time-differencebetween the two images and a warping function used to align the images.There are two main advantages of this approach over traditionalcross-correlation or cross-equalization approaches: (1) stability of themethod and (2) generation of time-continuous time-lags and warpingfunctions, a direct indicator of velocity change in the reservoir.However, once the S₁ and S₂ data are generated by the present inventivemethod, the cumulative S-wave time difference may alternatively becalculated by cross-correlation, cross-equalization, or any similarcomparison method, but using the similarity attribute is preferredbecause of, at least, stability and robustness.

One embodiment of the present disclosure uses a similar approach toregister S₁ and S₂ images and to produce a time-continuous S-wavetime-difference. As appreciated by those of ordinary skill, Fomelteaches an image registration technique to remove artifact time-lapsedifferences caused by velocity changes. Where Fomel aligns P-wave dataat different times, certain embodiments of the present disclosure alignsimultaneous S-wave data after it has been separated into its S₁ and S₂modes. Because of the time-continuous nature of the time-difference, theinventive methodology allows for the computation of the S-wavetime-difference between any two reflection events. Such information canlead surveyors to highly fractured zone in the subsurface.

In order to perform a similarity-based S-wave time-differencecalculation, the window-based S-wave time-difference applied in step 115on the registered slow S-wave volume is removed in order to generate anuncompensated slow S-wave (S₂) volume (117). A local similarity scanbetween S₁ and uncompensated S₂ volumes is then performed (119) to findthe optimum warping function, which is the name given to a function oftime that compensates the S₂ volume so that it is in phase with thecorresponding S₁ volume. The warping function also indicates the extentto which the cumulative S₂ time shift departs from a linear function oftime for a particular fracture system. The correct warping function (orwarping function trend), w(t), is generated by detecting the areas ofstrong similarity (121) and is calculated to minimize the differencebetween the S₁ and uncompensated S₂ volumes when applied to theuncompensated S₂ volume. In some embodiments, a smoothness constrain isapplied on the warping function while it is being calculated. Any kindof smoothness constrain me be applied on the warping function in time(or vertical) and/or lateral (X and Y) directions. In one embodiment,the user is prompted to provide the number of points in each directionto use in order to smooth out the warping function.

FIG. 2A is an illustration depicting a similarity scan (201) as thebackground and the chosen warping function (203) for a seismic tracepair in S₁ and uncompensated S₂ seismic sections.

Referring again to FIG. 1, based upon the generated warping function, acumulative S-wave time-difference is then calculated (123). In oneembodiment, the cumulative S-wave time-difference, ΔT_(S) ^(cum)(t), isdetermined using the following equation, which follows from thedefinition of the warping function:ΔT _(S) ^(cum)(t)=(w(t)−1)t,  (1)where t is the reflection travel time. FIG. 2B illustrates thecalculated cumulative S-wave time-difference (205) as a function of timebased on the warping function depicted in FIG. 2A.

Based upon the cumulative S-wave time-difference, an instantaneousS-wave time difference is calculated (125). The fracture density of asubsurface formation is directly proportional to the instantaneousS-wave time difference. In one embodiment, the instantaneous S-wavetime-difference, δt_(S) ^(ins)(t), is determined at each seismic datatime sample using the following equation:

$\begin{matrix}{{\delta\;{t_{S}^{ins}(t)}} = \frac{d\;\Delta\;{T_{S}^{cum}(t)}}{d\; t}} & (2)\end{matrix}$

The described process may then be performed for every trace pair in thevolume, thereby resulting in a cumulative S-wave time-difference volumeand an instantaneous S-wave time-difference volume being generated. Inanother aspect of the present disclosure, a new compensated S₂ volumemay be generated by using ΔT_(S) ^(cum)(t) volume to time-shiftuncompensated S₂ volume (127).

It is important to note that the steps depicted in FIG. 1 are providedfor illustrative purposes only and a particular step may not be requiredto perform the inventive methodology. As a non-limiting example, asurvey may decide not to calculate an instantaneous S-wave timedifference and also forego the generation of a new compensated S₂ volumebased upon survey objectives and/or computer processing speed or memory.Again, the claims, and only the claims, define the inventive system andmethodology.

EXAMPLE

To demonstrate applicability of the disclosed methodology, the presentinventors applied the method on synthetic seismic modeling data. FIG. 3shows S₁ traces (301) in black and uncompensated S₂ traces (303) in red.FIG. 3 illustrates that events in S₁ and uncompensated S₂ volumes arenot aligned in time.

FIG. 4 illustrates the S₁ traces (401) in black and the new compensatedS₂ traces (403) in red. As shown, the events in S₁ and compensated S₂volumes are now aligned. FIGS. 5 and 6 depict the Δ_(S) ^(cum)(t) andδt_(S) ^(ins)(t) volumes, respectively. This information allows anoperator or surveyor to compute an S-wave time lag and/or instantaneousS-wave time difference between any two reflection events. This new,additional information will help to locate highly anisotropic orfractured zones: critical information for well planning inunconventional reservoirs.

A review of FIGS. 5 and 6 shows that the fracture intensity of thesubsurface formation appears to be large in an area corresponding to atravel time between 3 and 3.5 seconds. More specifically, the area ofgreatest fracture intensity is identified at (501) and (601) in FIGS. 5and 6, respectively.

While the invention has been illustrated and described in detail in thedrawings and foregoing description, the same is to be considered asillustrative and not restrictive in character, it being understood thatonly the preferred embodiments have been shown and described and thatall changes and modifications that come within the wording of the claimsare desired to be protected. It is also contemplated that structures andfeatures embodied in the present examples can be altered, rearranged,substituted, deleted, duplicated, combined, or added to each other. Thearticles “the”, “a” and “an” are not necessarily limited to mean onlyone, but rather are inclusive and open ended so as to include,optionally, multiple such elements.

REFERENCES

-   1. Aguilera, R., “Naturally fracture reservoirs”, PennWell Book,    Tulsa, 1995.-   2. Nelson, R. A., “Geologic analysis of naturally fractured    reservoirs”, Gulf Publishing Company, Houston, 2001.-   3. Sinha, B. K., Norris, A. N. and Chang, S., “Borehole flexural    modes in anisotropic formations”, Geophysics 59, 1037-1052 (1994)-   4. Schoenberg, M. and Douma, J., “Elastic wave propagation in media    with parallel fractures and aligned cracks”, Geophysical Prospecting    36, 571-590 (1988).-   5. Ruger, A. and Tsvankin, I., “Using AVO for fracture detection:    Analytic basis and practical solutions”, The Leading Edge 16, 1429    (1997).-   6. MacBeth, C. and Crampin, S., 1991, “Comparison of signal    processing techniques for estimating the effects of anisotropy”,    Geophysical Prospecting 39, 357-386.-   7. Alford, R. M., “Multisource multireceiver method and system for    geophysical exploration”, U.S. Pat. No. 5,343,441 (1994).-   8. Winterstein, D. F. and Meadows, M. A., “Shear-wave polarization    and subsurface stress directions at Lost Hills field”, Geophysics    56, 1331-1348 (1991).-   9. Gaiser, J. E., “3-D converted shear wave rotation with layer    stripping”, U.S. Pat. No. 5,610,875 (1997).-   10. Tsvankin, I., “Seismic signatures and analysis of reflection    data in anisotropic media”, Pergamon, New York (2001).-   11. Thomsen, L., Tsvankin, I. and Mueller, M. C., “Coarse-layer    stripping of vertically variable azimuthal anisotropy from    shear-wave data”, Geophysics, 64, 1126-1138 (1999).-   12. Granger, P. Y., Bonnot, J. M., Gresillaud, A. and Rollet, A.,    “C-wave resolution enhancement through birefringence compensation at    the Valhall field”, Society of Exploration Geophysicist Annual    Conference, 2001.-   13. Crampin, S., “Evaluation of anisotropy by shear-wave splitting”,    Geophysics, 50, 142-152 (1985).-   14. Gaiser, J. E., “Application for vector coordinate systems of 3-D    converted-wave data”, The Leading Edge, 18, 1290-1300 (1999).-   15. Alford, R. M., “Shear data in the presence of azimuthal    anisotropy: Dilley, Tex.”, SEG Expanded Abstracts, 5, 476-479    (1986).-   16. Thomsen, L., “Converted-wave reflection seismology over    inhomogeneous media”, Geophysics, 64, 678-690 (1999).-   17. Haacke, R. R., Westbrook, G. K. and Peacock, S., “Layer    stripping of shear-wave splitting in marine PS waves”, Geophysical    Journal International, 176, 782-804 (2009).-   18. Bansal, R., Matheney M. and Liu, E., “True—amplitude    Layer-stripping in fractured media”, U.S. Patent Application No.    61/484,949.-   19. Fomel, S. and Jin, L., 2009, “Time-lapse image registration    using the local similarity attribute”, Geophysics, 74, A7-A11.-   20. Bakulin, A., Grechka, V., and Tsvankin, I., “Estimation of    fracture parameters from reflection seismic data—Part I: HTI model    due to a single fracture set”, Geophysics 65, 1788 (2000)

What is claimed is:
 1. A method, comprising: obtaining S-wave seismicdata of a subsurface formation, wherein the S-wave seismic data isacquired from the subsurface formation using multi-component receiversadapted to measure a plurality of motion vector components including twohorizontal components, the S-wave seismic data comprising a plurality ofseismic time samples; separating the S-wave seismic data into fast(“S₁”) and slow (“S₂”) data; using a computer to compute localsimilarity of the S₁ and S₂ data at each time sample; using a computerto compute at each time sample a cumulative time-difference by which theS₂ data lags the S₁ data from the local similarity, wherein computingthe cumulative time-difference includes, performing a local similarityattribute scan between the S₁ and S₂ data; generating a warpingfunction, w(t), by detecting areas of similarity in the local similarityattribute scan, and determining the cumulative time-difference, ΔT_(S)^(cum)(t), based upon w(t); calculating an instantaneous S-wavetime-difference, δt_(S) ^(ins)(t), at each seismic time sample from thecumulative S-wave time difference; estimating fracture density withinthe subsurface formation from the cumulative time-difference, ΔT_(S)^(cum)(t); and drilling a well into the subsurface formation based atleast in part on the fracture density.
 2. The method of claim 1, whereinthe cumulative time-difference, ΔT_(S) ^(cum)(t), is determined using aformula that can be expressed asΔT _(S) ^(cum)(t)=(w(t)−1)t, where t is reflection travel time.
 3. Themethod of claim 1, further comprising: generating time-compensated S₂data based upon the →T_(S) ^(cum)(t) to time shift the S₂ data.
 4. Themethod of claim 1, wherein the subsurface formation has a primaryfracture direction, the S₁ and S₂ data is generated by vector-rotatingthe horizontal components to directions parallel and normal to theprimary fracture direction.
 5. The method of claim 1, wherein the S₁ andS₂ data are generated by performing a layer stripping process on theS-wave seismic data to generate the S₁ data and a registered S₂ datawhich has a time-difference applied and removing the time-differenceapplied to the registered S₂ data to generate the S₂ data.
 6. The methodof claim 1, wherein the instantaneous S-wave time-difference iscalculated using a formula that can be expressed as${\delta\;{t_{S}^{ins}(t)}} = {\frac{d\;\Delta\;{T_{S}^{cum}(t)}}{d\; t}.}$7. The method of claim 6, further comprising: estimating the fracturedensity based on δ_(S) ^(ins)(t).
 8. The method of claim 1, furthercomprising: producing hydrocarbons from the well.